/*
 * Tema 2 ASC
 * 2020 Spring
 */
#include "utils.h"

/*
 * Add your optimized implementation here
 */


double* my_solver(int N, double *A, double *B)
 {
 	int i, j, k;
 	double *BAt, *AB, *C;
	register double *ptr1, *ptr2, *ptr3, *ptr4;
 	register double sum;
    register int index;

 	BAt = (double *) malloc(N * N * sizeof(double));
 	AB  = (double *) malloc(N * N * sizeof(double));
 	C   = (double *) malloc(N * N * sizeof(double));

 	for (i = 0; i < N; i ++) {
        index = i * N;
		ptr3 = &(BAt[index]);
		ptr4 = &(AB[index]);

 		for (j = 0; j < N; j++) {
 			sum = 0.0;
			ptr1 = &A[(j * N) + j];
			ptr2 = &B[(index) + j];

            /* B * (A transpus) */
 			for (k = j; k < N; k++) {
                /* sum += A[(j * N) + k] * B[(i * N) + k]; */
				sum += *ptr1 * *ptr2;

				ptr1++;
				ptr2++;
			}

			*ptr3 = sum;
 			/* BAt[(i * N) + j] = sum; */

 			sum = 0.0;
			ptr1 = &A[(index) + i];
			ptr2 = &B[(index) + j];

            /* A * B */
 			for (k = i; k < N; k++) {
				sum += *ptr1 * *ptr2;
 				/* sum += A[(i * N) + k] * B[(k * N) + j]; */
				ptr1++;
				ptr2 += N;
			}

            /* AB[(i * N) + j] = sum; */
			*ptr4 = sum;

			ptr3++;
			ptr4++;
 		}
 	}

 	for (i = 0; i < N; i++) {
        index = i * N;
        ptr3 = &(C[index]);
        ptr4 = &(BAt[index]);

 		for (j = 0; j < N; j++) {
 			sum = 0.0;
			ptr1 = &(A[(i * N) + i]);
			ptr2 = &(AB[(i * N) + j]);

 			for (k = i; k < N; k++) {
                /* sum += A[(i * N) + k] * AB[(k * N) + j]; */
				sum += *ptr1 * *ptr2;

				ptr1++;
				ptr2 += N;
			}

            /* C[(i * N) + j] = BAt[(i * N) + j] + sum; */
            *ptr3 = *ptr4 + sum;

            ptr3++;
            ptr4++;
 		}
 	}

 	free(BAt);
 	free(AB);

 	return C;
 }
